#args <- commandArgs(TRUE)
library(MASS)
library(sfsmisc)
snp<-read.table("/data1/bsi/BORA_processing/devel/eqtl/parallelize_genotypingsurvival/Matrix_eQTL_R_Source/SNP.txt",sep="\t",head=T)
gene<-read.table("/data1/bsi/BORA_processing/devel/eqtl/parallelize_genotypingsurvival/Matrix_eQTL_R_Source/OGE.txt",sep="\t",head=T)
k1<-c("gene","snp","Intercept_Estimate","snp_Estimate","Intercept_Std.Error","snp_Std.Error","Intercept_t value","snp_t value","ROBUST_F_TEST_RLMFIT")
data2<-read.table("/data1/bsi/BORA_processing/devel/eqtl/parallelize_genotypingsurvival/Matrix_eQTL_R_Source/rlmlist.txt",sep="\t",head=T)
colnames(data2)<-c("rsid","geneid")
zz <- file("/data1/bsi/BORA_processing/devel/eqtl/parallelize_genotypingsurvival/Matrix_eQTL_R_Source/rlmlist_out.txt","w")
writeLines(paste(k1,collapse=" "),con=zz,sep="\n")
#sorting accordind to file order
#sort_order<-read.table(args[3],head=T)
#sort_order<-as.vector(t(sort_order))
#gene<-gene[,sort_order]
rownames(gene)<-gene$geneid
z<-1
while(z<nrow(data2))
{
	zz <- file("/data1/bsi/BORA_processing/devel/eqtl/parallelize_genotypingsurvival/Matrix_eQTL_R_Source/rlmlist_out.txt","a")
	z1<-z+100000-1
	if(z1 > nrow(data2))
	{
		z1<-nrow(data2)
	}
	data<-data2[z:z1,]
	colnames(data)<-c("rsid","geneid")
	data<-merge(data,snp,by="rsid")
	data<-merge(data,gene,by="geneid")
	i<-1
	while(i<nrow(data))
	{
        	print(i)
		#g2<-as.numeric(data[i,3:209])
        	#g1<-as.numeric(data[i,210:416])
		kk<-c()
		flm<-rlm(as.numeric(data[i,210:416]) ~ as.numeric(data[i,3:209]), family=gaussian(link="log"))
		#flm<-rlm(g1 ~ g2, family=gaussian(link="log"))
		kk<-c(as.character(data[i,2]),as.character(data[i,1]),as.vector(coefficients(summary(flm))),f.robftest(flm)$p.value)
		writeLines(paste(kk,collapse=" "),con=zz,sep="\n")
		i<-i+1
	}
	close(zz)
	z<-z+100000
}
